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Wave front pinning and propagation in damped chains of coupled oscillators are studied. There 
are two important thresholds for an applied constant stress F: for |F| < Fed (dynamic Peierls stress), 
wave fronts fail to propagate, for Fed < \F\ < Fes stable static and moving wave fronts coexist, and 
for |F| > Fc3 (static Peierls stress) there are only stable moving wave fronts. For piecewise linear 
models, extending an exact method of Atkinson and Cabrera's to chains with damped dynamics 
corroborates this description. For smooth nonlinearities, an approximate analytical description is 
found by means of the active point theory. Generically for small or zero damping, stable wave front 
profiles are non-monotone and become wavy (oscillatory) in one of their tails. 

PACS numbers: 45.05.+x; 05.45.-a; 83.60.Uv 



I. INTRODUCTION 

Wave fronts and pulses play important roles in many 
physical systems. Examples abound: the motion of dislo- 
cations J ^ or cracks ^3] in crystalline materials, atoms 
adsorbed on a periodic substrate 4] , the motion of elec- 
tric field domains and domain walls in semiconductor su- 



pulse propagation through myelinated 



perlattices 

nerves j?!, pulse propagation through cardiac cells 
etc. Furthermore, these localized waves often play an 
important role in Statistical Mechanics or Quantum 



critical value 11 ll . Wave front pinning in discrete sys- 
tems may be related to such different physical phenom- 
ena as the existence of Peierls stresses in continuum me- 
chanics ^12! or the relocation of electric field domains in 
semiconductor super lattices In the continuum limit, 
the width of the pinning interval (range of the exter- 
nal force for which wave fronts fail to propagate) tends 
to zero exponentially fast and many authors have calcu- 
lated the critical force for different models in this limit 
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Field Theory When wave fronts or pulses are solu- 
tions of spatially discrete systems, they often fail to prop- 
agate unless an external force or parameter surpasses a 



Not surprisingly, wave front motion and pinning are 
different depending on the dynamics describing the 
model at hand. To be precise, let us consider a chain 



of nonlinear oscillators, diffusively coupled and subject 
to an external force F which acts as a control parameter: 

m— -5- + -3— = Un+l - 2Un + Un-1 - A g{Un) + F. (1) 

Typical nonlinearities g{u) = V'{u) are cubic, such that 
Ag{u) - F has three zeros, Ui{F/A) < U2{F/A) < 
U^{F/A) in a certain force interval {g'{Ui{F/A)) > 
for z = 1,3, g'{U2iF/A)) < 0). Moreover, g{u) is sym- 
metric with respect to U-){0). Examples are the over- 

n 

damped Frenkel-Kontorova (FK) model {g = sinu) l] 
and the quartic double well potential {V = {u^ — 1)^/4) 
\l\ . A > measures the strength of the coupling and m 
the relative strength of inertial and friction terms. Wave 
front solutions m„ = w{n — ct) join the two stable con- 
stant states Ui{F/A) and U3{F/A) (or viceversa) as n 
increases from —00 to 00. 

Consider the extreme cases of conservative (m — 00) 
and overdamped dynamics (m = 0). In the overdamped 
case, wave fronts generically either move if |_F| > _Fc > 
or are pinned if \F\ < Fc The dcpinning transition 
at Fc was described by Carpio and Bonilla [11 1 for large 
and moderate values of A, by King and Chapman 15] in 
the continuum limit A — > 0, and by Path 16] for a piece- 
wise linear g{u). In the conservative case and for generic 
cubic nonlinearities g{u), there are two critical forces Fed 
and Fes with < Fed < Fes — Fe- Wave fronts may 
propagate stably for \F\ > Fed but there are stable sta- 
tionary (pinned) wave fronts if \F\ < Fe- Thus pinned 
and moving wave fronts may coexist if Fed < 1^1 < Fes- 
The values Fes and Fed correspond to the static and dy- 



namic Peierls stresses of the literature on dislocations 
Atkinson and Cabrera found exact expressions for the 
wave fronts corresponding to a piecewise linear g{u) and 
calculated the relationship between F and wave front ve- 

n 

locity |17| . An approximate theory was found somewhat 
earlier by Weiner |l^ . More recently, Schmidt ^ and 
later authors [2^ L J ] found exact monotone wave fronts 
of conservative or overdamped systems by constructing 
models with nonlinearities such that the desired wave 
fronts were solutions of the models. In particular, Flach 

n 

et al |21| showed coexistence between moving and pinned 
monotone wave fronts of a discrete system with conser- 
vative dynamics and _F = 0. On the other hand, for a 
sine nonlinearity and F — 0, the numerical computations 
of Peyrard and Kruskal 22] show that an initial profile 
close to the continuum sine Gordon soliton loses energy 
via emission of phonons and it becomes pinned after a 
sufficiently long time interval. They also found stable 
moving wave fronts for small positive F, consistent with 
our previous statement that, generically. Fed > 0. 

In this paper, we study the wave fronts of the damped 
system and the transitions at Fed and Fes- In con- 
trast with previous work, we find wave fronts that are 
non-monotone, presenting wavy tails at one or both sides 
of a transition region in which the profile w{n—cT) jumps 
an amount close to [UsiF/A) - Ui{F/A)]. We call them 
wavy wave fronts. These fronts with wavy profiles per- 
sist even in the conservative limit (m 00) 2^, and in 
fact Atkinson and Cabrera's wave fronts are also wavy. 
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as these authors would have found out had they depicted 
their exact expression graphically. In the overdamped 
limit m — > 0, F^d Fc and the wave front profiles be- 
come monotone. We have thus arrived to a general pic- 
ture of wave fronts in discrete chains of coupled nonlinear 
oscillators with m > 0. 



The rest of the paper is organized as follows. Section 
im considers Eq. Q| with a piecewise linear g{u). We 
find exact formulas for the wave front profiles in the gen- 
eral damped case following the method of Atkinson and 

n 

Cabrera's [17| . These profiles are often wavy and they are 
asymptotically stable in the damped case. It is important 
to obtain them for two reasons: (i) there are very few ex- 
act wave front solutions that are non-monotone, and (ii) 
in the limit of large inertia, it is hard to discriminate 
numerically between wavy wave fronts traveling with dif- 
ferent velocities or having different profiles. Exact so- 
lutions make good benchmarks for numerical methods. 
The results for the damped model with a generic cubic 
nonlinearity are presented in Section IIIII We calculate 
the static and dynamic Peierls stresses for typical values 
of A and m. A characterization of these stresses is given 
in terms of our active point theory. Sec. IIVI contains a 
discussion of our results. 



II. EXPLICIT CONSTRUCTION OF WAVE 
FRONT PROFILES 

Let us rescale time in Eq. so that t — t j ^/m^ and 
consider a piecewise linear g(u): 
d Un du 

+ = Un+l - 2u„ + Un-l - Ag{un) + F, (2) 

{Un + l, for Un < 0, 
(3) 
Un - 1, for Un > 0, 

where a = l/y/m. Notice that g(u) = u + 1 — 2H{u), 
where H{x) = 1 for x > and H{x) ~ for x < is the 
Heaviside unit step function. Let us consider a smooth 
wave front profile w„ — w{x) = v{x) ~ 1, x = n — ct, 
moving rigidly with velocity c. We center the wave front 
so as to have w{0) — 0. Taking into account that g{u) 
is an odd function and using the front profile = 
w{n — ct), we can see that the following transformations 
leave Eq. (jSJ for w{x) invariant: 

{x, w, c, F) {—X, w, — c, F) 
(x, w, c, F) {x, —w, c, ~F) 
{x, w, c, F) — > {—X, —w, — c, —F). (4) 

Let us consider now the case F > 0, c > and ^'(O) < 0, 
i.e., a wave front profile that decreases in the transition 
region about a: = 0. The transformations Q yield a 
profile with: (i) c < 0, -F > increasing in the transition 
region, (ii) c > 0, F < increasing in the transition 
region, and (iii) c < 0, F < decreasing in the transition 
region. Thus we find that sgnw — — sgn(xcF), g{w) — 
w + 1 — 2 xsgn(cF)), and we can restrict ourselves 
to considering the case F > 0, c > and w'{0) < 0: 
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FIG. 1: Symmetries in the wavefront solutions for A = 0.25, 
a = 0, c = 0.5. and F = 0.009. 

All other three possible cases can be obtained from our 
results by using Eq. Q); see Figure 

The wave front profile v{x) = w{x) + I satisfies: 

c^v"{x) - acv'{x) - [v(x + f ) - 2v{x) + v{x - 1)] 

+Av{x) = 2AH{-sgn{cF)x) + F, (5) 

with i;(0) = f . We can calculate v{x) by using the con- 
tour integral expression for the step function: 



I 



2711 Jq k 



■dk. 



(6) 



Here C runs over the real axis in the complex k plane 
passing above the pole at fc = as in Fig. |21 For x > 
(resp. a; < 0), C is closed by a semicircle in the upper 
(resp. lower) half plane oriented counterclockwise (resp. 
clockwise). 

Then the solution of Eq. © is 

, , F A [ expfifc sgn(ci^)a;] dk 

"(a^) = T - Z7 / rm—^ ' (7) 




Re(k) 



A TTi J (J kL{k,a) 



FIG. 2: Contour for the Heaviside step function @ and the 
integral formula 10 -© when a ^ 0. 



L{k,a) = A + 4sin2 ( ^ ) - k^c^ - ik\c\a sgn{F) . (8) 



All the zeros of the function L{k, a) given by Q are 
complex for a > 0, and they correspond to exponentially 
localized modes. The nonzero poles of the integrand in 
Eq. can be found graphically by plotting the curves 
Rei(fc, a) = and ImL(fc, a) = in the complex fc-plane, 
as depicted in FigO When a 0, a finite number of 
poles tend to the real axis, whereas infinitely many keep 
a nonzero imaginary part even at a = 0. The poles on 
the real axis correspond to radiation modes, cause os- 
cillations in the wave front tails, and their number in- 
creases as c decreases; see Figs. EJa) and (b). The purely 
imaginary poles of Figures O and ^c) yield the central 
monotone part of the wave front profiles. For a — 0, the 
integration contour in Eq. Q avoids poles on the real 
axis according to a criterion due to Atkinson and Cabrera 

n . . 

[17 1, shown in Figure |S1 and derived later on this Section. 
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FIG. 3: Complex poles as intersection of the curves 
ReL(fc, a) = and ImL(fc, a) = when c = 0.2, A = 0.25 
and: (a) q = 0, (b) a = 0.1, (c) a = 1. 

We will use {T)) and the method of residues to construct 
profiles satisfying v{x) > 1 for a; < and v{x) < 1 for 
X > 0. Notice that we can obtain a complex dispersion 
relation between to — kc and k from L{k,a) = 0. The 
contour choice and the fact that a > give rise to an ex- 
ponential decay of v{x) to its asigned values at a; = ±oo. 
When a — 0, the wave fronts may exhibit undamped 
oscillations extending all the way to infinity. 

The condition v{0) = 1 yields a relationship between 
the wave front velocity c and the external force F: 



F A r dk 

A Tri J (J kL{k, a) 

^- y 

L{p,a}=QJm{p)>a^ ' 



2A 



(9) 



where we have assumed that cF > and v'{l) < 0. The 
resulting function F(c) can be calculated by computing 
this series of residues numerically. Once F{c) is known, 
Eq. ^ can be used to compute the wave front profiles for 



FIG. 4: Real poles in the conservative case q = with A = 
0.25: (a) c — 0.5 (b) c = 0.01; Purely imaginary poles when 
A = 0.25 (c) c = 0.5, a = 1 
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FIG. 5: Contour for the integral formula Q-© in the con- 
servative case Q = for c > and F > 0. 



a pair (c, F(c)). We shall now show how this construction 
works out for a = 0, a = oo and for a finite. 



A. Conservative case: a = 

It is instructive to see what happens in the conser- 
vative hmit a 0+. In this case, L{k,Q) = has 
real solutions, and we need a criterion to move the con- 
tour C above or below the corresponding poles in the 
integral Q. To obtain it, we shall use the notation 
L{k) — L{k,0). Let ko be a real zero of L{k). The 
complex zero of L{k, a) that becomes fco at a = satis- 
fies = Lk{ko,0) {k — kQ) + La{ko,0) a+ . . ., which yields 
k ^ kg — La{kQ, Q)a/ Lk{kQ, 0), that is. 



k ^ kn 



iafco|c|sgni^ 



L'iko) 

sgn(Imfc) = sgn(cF) sgn[ckoL' (ko)], (10) 

as a — > 0-I-. We know that the contour C in Eq. © lies 
in the upper (resp. lower) half plane provided cFx > 
(resp. cFx < 0). Therefore, the poles whose residues 
count must satisfy cFxlmk > 0. Then Eq. H1U|) implies 
that we should count poles satisfying xckQL'{ko) > 0. 
The physical meaning of this criterion becomes clear if 
we calculate the group velocity corresponding to mode 
ko. L{k) = yields (fcc)^ = [uj{k)]^ = .4 + 4 sin2(fc/2). 
Then 2w(fc)cj'(fc) = 4sin(fc/2) cos(fc/2) = L'{k) + 2kc^. 



Thus Vg = Lu'{k) obeys 



L'jk) 
2kc 



sgn(wg - c) = sgn[cA:L'(fc)], (11) 



and we observe that the poles whose residues contribute 
to the solution satisfy {vg ~ c)x > 0. This was the crite- 
rion used by Atkinson and Cabrera [ij: All modes with 
Vg > c must appear ahead of the wave front (x > 0), 



all those with Vg < c must appear behind (x < 0). See 
Figure El 

For c > and a = 0, the condition ^ becomes |l7 

2A 



F 
A 



E 



k\L'{k)\- 



(12) 



L(fc)=0,fc>0 

This formula follows straightforwardly from the fact that 
L(fc) and kL'ik) are even functions of (real) k and sym- 
metry considerations. Notice that our assumption c > 
has yielded F > 0. The relation F{c) given by (|12|l is 
plotted in Figure El for a value A — 0.25 (see also Fig. 
3 of Ref. Q). For a given value of the external force 
F, there may be several values of admissible velocities c, 
each corresponding to a different wave front profile. Thus 
different families of wave fronts (not all of them stable) 
may coexist for the same value of F. The function F{c) 
presents different vertical asymptotes at positive values 
Ci, ci > C2 > . . ., where both L{k) and L'{k) and van- 
ish for positive k. We have c„ ~ vG4/(27rn) as n ^ oo, 
so that the vertical asymptotes accumulate at c = as 
suggested by Fig. 01 In fact, the velocity can be elimi- 
nated from the two conditions L{k) = L'{k) = yielding 
yl -I- 2 (1 — cosfc) — ks'mk. For large values of k, this 
gives k 27171 + A/{2T:n), as the integer n — s- oo. The 
condition L'{k) = then yields the previous formula for 

The range of physical interest corresponds to ci < 
c < 1 (wave front velocities larger than the largest res- 

n 

onant velocity but smaller than the sound speed |17j). 
Let {cjn, F(c„i)) be the minimum of the first and fastest 
branch of F{c). Then Fed = F{crn) yields the dynamic 




FIG. 6: (a) F{c) for a = 0, A = 0.25 as computed from 
formula IIH , (b) Zoom in the region of resonances for c small. 
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FIG. 7: Wave front profiles for A = 0.25, a = and: (a) 
c = 0.23, (b) c = 0.2. 

Peierls force, under which the 'physical' branch of stable 
wave front profiles ceases to exist. This force is smaller 
than the static Peierls force, — Ai (A + 4)^i, for the 
piecewise linear model |l7|. 

Using the radiation condition [vg — c)x > 0, we have 
plotted in Figuresnand[7|several wave front profiles. Fig- 



ure a) shows a wave front profile for a = 0, A = 0.25 
and F = 0.009. The wave velocity is c = 0.5. The pro- 
file has been numerically approximated by computing the 
contour integral Q as the series of residues truncated to 
a few terms. It is interesting to observe that the right 
tail of the wave front decays fast to Ui{F/A) whereas 
the left tail oscillates about Uj,{F/A) with uniform am- 
plitude. The reason for this behavior is that for the above 
parameter values, there are two real zeros of L(k) at zt/co 
with k^L' [kf)) > that contribute a non-decaying oscil- 
lation to the left tail. An infinitesimal amount of friction 
would dampen these oscillations by contributing a multi- 
plicative factor exjp{~[akQc/ L' (ko)]} to their amplitude. 
Figure [T^a) shows a wave front profile for c = 0.23. No- 
tice the decaying small oscillation in the right tail. We 
have still two real poles, but now c — 0.23 is placed at 
the left of the minimum in the first branch of F(c); see 
Figure Ela). Figure [7{b) shows a wave front profile for 
c = 0.2, in the second branch of F{c), past the first res- 
onance. Now we have three pairs of real poles. Two 
of them contribute to the oscillation in the left tail, the 
other one produces the oscillation in the right tail. 

The dynamical stability of the constructed wave front 
solutions can be numerically checked by using their com- 
puted profile as initial data to solve Eq. ^ with a = 

1 

24| . The results are compared at a fixed time t = 60 to 
the expected configuration w{n — ct) in Figure |S1 The 
choice c = 0.5 seems to produce a stable wave front. The 
choice c = 0.23, (still on the first and fastest branch of 
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FIG. 8: Dynamical stability when A = 0.25, a = 0. We 
compare u„{t) (solid line) to w{n — ct) (dot-dashed line) for 
(a) c = 0.5, (b) c = 0.23, (c) c = 0.2. 



FIG. 9: Overdamped limit when A — 0.25. (a) Wave front 
profile for c = 0.02. (b) F{c) as computed from formula @. 

B. Overdamped limit: a = cxd 



F(c) but to the left of its minimum), evolves towards a 
static front. The choice c ~ 0.2 (on the second branch of 
F(c)), evolves towards a wave front moving faster than 
expected, with a speed on the first branch of F(c). Thus 
our numerical results seem to indicate that stable wave 
fronts have velocities on the first and fastest branch of 
F{c) with F'{c) > 0, to the left of the minimum speed 
on this branch, c„i > 0. Then Fd — F[cm), and stable 
wave fronts with v'(\) < have speeds larger or equal 
than Cm- 

We have found wave front profiles with oscillatory tails 
that seem stable under small disturbances. One question 
that comes to mind is whether these profiles occur in 
models with smooth nonlinearities. The answer is yes: 
See an explicit construction in Appendix A. 



The results in the overdarnped limit m = are consis- 
tent with previous work llSl Il6| : there are one wave 
front profile and one c for each fixed F above a threshold, 
Fc- Wave front profiles are monotone, and they resemble 
staircases for c small. See Figure |S| 

C. Finite damping: a > 

The results for finite damping interpolate between the 
conservative and overdamped cases. For small a, the 
function F{c) and the wave front profiles are non mono- 
tone although their oscillations decay as n ^ ±oo; see 
Fig. ^1 Fig. im shows a comparison between Un{t) for 
a = (for the same values as in Fig. ISJ. We observe that, 
for this small damping, the corresponding wave fronts 
have the same stability properties as in the conservative 
case: dynamically stable for c > Cm, and unstable for 
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FIG. 10: Results for A = 0.25, a = 0.01 (a) Wave front profile 
for c — 0.5, (b) -F(c) as computed from formula @. 
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FIG. 11: Comparison between calculated for a = for 

the same values as in Fig. |Hl(solid line) and u„{t) for a = 0.01 
(dot-dashed line), (a) F = 0.009, (b) F = 0.01, (c) F = 0.029. 

c < Cm- Moreover there are dynamic and static Peierls 
stresses which are different from each other, as in the case 
a = 0. 

It is interesting to ascertain the shape of the function 
F{c) because its inverse function constitutes a bifurcation 
diagram clarifying wave front depinning. In the conser- 



vative case, F{c) has infinitely many vertical asymptotes 
(resonances) accumulating at c = 0. As a ^ 0+, the 
conditions L{k, a) ~ Lk{k, a) = yield 

akoc 



k 
X 



ko + 



cos fco — 
sgn[(cos fco — c^)fcoF] 



e'^^+'''\ 1 = 0,1, (13) 
(14) 



This formula shows how the double poles fc = fco of the 
conservative case split when an infinitesimal friction is 
present. Then the vertical asymptotes of F(c) at c = c„ 
give rise to local maxima of F{c) for small a > 0. These 
maxima are hard to resolve numerically (see Fig. Iiur b)). 
but they can be approximately calculated as follows. 
Only poles with positive imaginary part contribute to 
the sum in Eq. 0. For these poles, the exponential fac- 
tor in Eq. H13|) is {x + i)/V2, and their contribution to 
the sum in Eq. is approximately given by 

fco sgn(cosfco — c^) |fcoc(cosfco — c^)|2 
To this expression, we should add its complex conjugate, 
a contribution to the sum in Eq. (|5J due to the pole —fco. 
If we keep only these contributions in Eq. (jS)), thereby 
assuming that the considered maximum of F{c) is large, 
Fm 1, we obtain 



\F 



V2A^ 



M\ 



(15) 



|fco| 2 |c (cos fco — c^)! 2 
Now, in the conservative case, c„ ~ \/A/ (27rn) and fco ^ 
27m -I- A/{2TTn) as the integer n oo. Then the right 
side of Eq. H15|) becomes proportional to c„. When n 
is so large that c„ is no longer large compared to ^/a, 
other terms of the sum contribute appreciably to F in 



10 



the formula (jS)). We conjecture that these contributfons 
add to Fes, 



\Fm\~Fc: 



^ ^5 c„, c„ ^ ^, (16) 



so that the maxima of F{c) accumulate near c = as the 
integer n oo. We have depicted schematically the re- 
sulting F{c) and the bifurcation diagram of c versus F in 
Fig. El Corresponding to the infinitely many local ex- 
trema in Fig ll2r a). there are infinitely many limit points 
(saddle-node bifurcations) in Fig. I12r b'l. Our numeri- 
cal results indicate that only the branch of wave fronts 
with larger velocities in the physical interval Cm < c and 
1^1 > Fed are stable. This can be understood from the 
factorization theorems in According to these the- 

orems, one eigenvalue of the linear stability problem cor- 
responding to the solution branches in the bifurcation 
diagram changes sign at limit points. If we use that the 
branch of wave fronts with larger velocity is stable, the 
branch that coalesces with it at the limit point F — Fed 
is unstable. This branch coalesces with another one at 
another limit point with larger F, and, there, a different 
eigenvalue of the linear stability problem changes sign 
from negative to positive. If this is so, the new branch is 
also unstable and all other wave front solution branches 
in Fig. I12r b) could also be unstable. 

For larger values of a, the wave front profiles become 
monotone, the oscillation amplitudes in F{c) decrease 
and become difficult to appreciate; see Fig. E| The tran- 
sition from one parameter range to the other one occurs 
when the contribution from poles with small imaginary 
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FIG. 12: (a) Schematic function F{c) for small a > showing 
infinitely many maxima accumulating at c = and F = Fes ■ 
(b) The bifurcation diagram of wave front velocity versus F: 
there are infinitely many limit points (saddle-node bifurca- 
tions) corresponding to the extrema of F{c) in the interval 

Fed < F < Fes- 
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FIG. 13: Results for A = 0.25, q = 1. (a) Wave front profile 
for c — 0.5, (b) F{c) as computed from formula Q. 



part in Eq. becomes relevant. 




FIG. 14: Results for the Frenkel-Kontorova model with A = 
10: (a) Dynamic (solid line) versus static (dot-dashed line) 
thresholds as functions of a; (b) numerical velocities as func- 
tions of F for decreasing values of a: 3 (circles), 1.5 (squares), 
1 (asterisks), 0.7 (triangles), 0.57 (crosses). 



FIG. 15: Same as Fig. ll4l for the model with a quartic poten- 
tial. In (b), the values of a are: 3 (circles), 1 (squares), 0.7 
(asterisks), 0.3 (triangles), 0.2 (crosses), 0.1 (diamonds). 

A. Numerical results 



III. WAVY WAVE FRONTS FOR GENERIC 
CUBIC NONLINEARITIES 



For generic smooth cubic nonlinearities g{u), we can- 
not construct the wave front profiles by using contour 
integrals. However, we can extend our previous theory 
of the active points ^ | for threshold phenomena to the 
case of finite damping. Thus we shall present numerical 
results for (relatively) large A showing that wave fronts 
are similar to those for piecewise linear g{u). Near Fc, 
we shall use the theory of active points to interpret nu- 
merical results. 



Figures FhI and IT^ show the dynamic and static Peierls 
stresses and the wave front velocity for the damped FK 
and cubic models, respectively. We have constructed the 
traveling wave fronts by solving numerically Eq. |5J in 
large lattices. For F > Fes, we choose as initial data 
a static step-like profile: it„(0) = U^^^F/A) if n > 0, 
u„(0) = U^^HF/A) if n < 0, u;,(0) = Vn. We use the 
stable zeros C/^^^ (F/A) and U^^'> (F/A) as boundary con- 
ditions for large \n\. The numerical solution Un{t) evolves 
very fast to a traveling wave = w{n — ct) with a 

fixed constant value for the speed c. For F below the 
static threshold, we choose as initial data the traveling 
solutions already found. As boundary condition, we use 
again U'-^HF/A) and U'-^^F/A). The numerical solu- 
tions Un{t) evolve to a traveling wave u„(i) — win — ct) 
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with a profile and speed adjusted to the new value of F, 
provided F is larger than the dynamical threshold Fd- 
Below that value, the waves are pinned. The behavior of 
c near Fd can be guessed from the known fact that the 
function F{c) has a minimum Fed on its fastest branch 
at c = Cm for piecewise linear g{u); see Fig. El Near this 
minimum, F ^ Fcd + J (c~c„iY with 7 > 0, as indicated 
in Fig.IlSta). This yields (c - Cm) - [{F - Fcd)ll]^, a 
scaling that can be seen in Figs. Iior b') and (c), corre- 
sponding to smooth g{u). The number 7 can be fitted 
by taking careful numerical measurements near F^d- This 
seems to provide a good fitting over an interval of stresses 
that increases as a decreases. For larger friction values, 
taking values of F farther from Fed produces a better fit 
to a scaling with the same exponent 1/2 but with differ- 
ent 7; see the squares in Figs. llSf bl and (c). 

Our numerical measurements of the speeds near Fed 
seem to indicate that (except in the overdamped limit 
771 = 0) there is a critical non zero speed c„i > below 
which front propagation cannot be sustained. In the co- 
existence region, Fed < F < Fes, shown in Figs. 1141 and 
1151 both the traveling wave fronts and the static wave 
fronts are dynamically stable. The wave front profiles for 
different damping values and the cubic g{u) are depicted 
in Fig. El These profiles oscillate more and more as the 
damping coefficient decreases. For sufficiently large a, 
the wave front profiles are monotone and become simi- 




Q 



lar to those calculated in the overdamped limit [11|. Fed 
and Fes are almost equal. Similarly, as A decreases and 



FIG. 16: (a) Function F{c) for the piecewise linear g{u) with 
A — 0.25 and a — 0. (b) Velocity versus applied stress for 
g{u) — sinu with A = 10 and a = 1.4 (squares) and a = 
0.57 (triangles), (c) Velocity versus applied stress for g{u) = 
u{u^ — 1) with j4 = 10 and a = 1.0 (squares) and a = 0.1 
(triangles). 

we approach the continuum limit, the gap between static 
and dynamical thresholds is difficult to appreciate. 

There is an important difference between models with 
a periodic nonlinearity such as FK and models with a cu- 
bic g{u). In both cases, wave fronts can be constructed 
numerically for sufficiently large values of the damping. 
For A = 10, wave fronts of the cubic model can be numer- 
ically found at least for a > 0.08, whereas those of the FK 
model are found easily for larger damping, a > 0.5. For 
smaller values of the damping, the amplitude of the wave 
front oscillatory tails becomes so large that the FK wave 
front profile tends to jump between different periods of 
the nonlinearity: moving staircases are thus generated. 
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FIG. 17: Wave front profiles for the quartic potential as seen 
from the trajectory of a single point u„{t): (a) a — 3, (b) 
a = 1.5, (c) a = 0.3, (d) a = 0.15. 

B. Active point theory 

To get approximate formulas for the wave front profile 
and velocity in the strongly discrete limit, A 3> 1, we 
can resort to the active point theory ^|. In this limit, 
there is one active point, say uo{t), and all others obey 
either u„ Ui{F/A) (for n > 0) or u„ UsiF/A) (for 
n < Q). We assume that the wave front we will construct 
has F > 0, c > and w'{0) < 0, as in the previous 
Section. According to Eq. the active point satisfies 
the approximate equation: 



duo ^ fF\^rr 



-2uo- Ag{uo) + F. 



(17) 



This equation has three stationary solutions for F < Fes, 
zi{Fc/A) < Z2{Fc/A) < Z3{Fc/A), (zi and are stable 
and Z2 is unstable), and only one stable stationary solu- 
tion, Z3, for F > Fes- The critical field Fc — Fes is such 



that the expansion of the right hand side of H17() about 
the two coalescing stationary solutions has zero linear 
term, 2 + Ag'{uo) — 0, and 



2uo + Ag{uo)^Ui[^)+U3(^]+Fe. (18) 



For F sligthly above Fc, wo(i) = uo{A, Fc) + vo{t) obeys 
the following equation: 



_ + a—^aiF-Fc) + bvo, 



a = 1 



1 



1 



Ag'iUi{Fc/A)) Ag'iU^iFc/A)) 



b = -^<?"(uo) >0, 



(19) 
> 0, (20) 
(21) 



where we have used 2 + Ag'{uo) — 0, (fll^ and ignored 
higher order terms. This equation has two distinguished 
limits, a < {F - Fc)i < 1 and a > {F - Fc)i . In 
the latter case, we can ignore the inertia in Eq. 119|l . 
The resulting reduced equation is exactly that analyzed 

n 

in Ref. 11 except for a trivial rescaling of the time. The 
solution of that equation blows up at times (t — io) ~ 



±TTa/[2y/ab {F - Fc)] {to is an arbitrary constant). Then 
the wave front velocity is approximately given by the 
reciprocal of the interval between two consecutive blow 
up times, namely, 

^ab (F - Fc) 



(22) 



After blow up, the wave front profile is reconstructed 
by inserting an inner layer, in which uo(t) obeys Eq. 
(|17|l with F = Fc, and it jumps from a neighborhood 
ofzi{Fc/A) to Z3{Fc/A) 11]. 

If a < (F-F,)4 < 1, we can ignore friction in Eq. 
(|19|l thereby obtaining a conservative dynamical system. 
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Eq. lUni) with a = 0, as our reduced equation. Its tra- 
jectories also blow up and the wave front velocity can be 
straightforwardly calculated as 



(23) 



where K{1/^) = 1.854075 is the complete elliptic inte- 
gral of the first kind with module k — l/\/2 ^261]. Con- 
trary to the overdamped case, a consistent inner layer 
connecting blowing up trajectories of the reduced equa- 
tion and trajectories of Eq. (|17|l for F = Fc and a = 
does not exist. This points out to a breakdown of the 
active point theory as a — > 0+, which is consistent with 
our conjectured bifurcation diagram in Fig. I12r b') for the 
piecewise linear model. Assuming that Fig. I12f b'l is also 
the bifurcation diagram of models with smooth nonlin- 
earities, a succession of infinitely many saddle-node bifur- 
cations (that accumulate at c = 0, F = Fes) connect the 
branch of stable wave fronts and that of stationary front 
solutions between Fed and Fes- If this is the case, seeking 
a description in terms of standard normal forms as H19() 
and scaling is rather hopeless. Not surprisingly, the sus- 
picious scaling (|23|) is hard to check numerically (the sta- 
ble branch of moving wave fronts ends at F = Fed < Fe), 
whereas the scaling (|22|) can be easily checked {Fed ~ Fc 



in the overdamped limit) 

Let us now come back to Eq. (|17l) to explain the 
coexistence of moving and stationary wave fronts for 
Fed < F < Fe. In Fig. ^ we have depicted the 
phase plane associated to Eq. (|17|l for the cubic non- 




FIG. 18: Phase plane of Eq. H17|l with the cubic nonlinearity 
for A = 10, a = 0.5 and F = 2.45. 
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FIG. 19: Comparison of Fcda = Fc2 and Fed as functions of 
the damping a for A = 10 and: (a) the cubic nonlinearity, (b) 
its relative error, (c) the FK nonlinearity, and (d) its relative 
error. 

linearity g{u) = u{v? — 1). The profiles in Fig. ITTI 
correspond to one active point jumping from a neigh- 
borhood of Ui{F/A) to a neighborhood of U:i{F/A). 
The right side of Eq. H17(l has three zeros in the inter- 
val [Ui{F/A),U:i{F/A)\ for F < F^, corresponding to 



two stable spiral points or centers (zi{F/A),0), i = 1,3, 
and one saddle point {z2{F/A),0). At F — Fc^ two 
of the zeros coalesce and for F > only the spi- 
ral point (z3(F/yl),0) remains. For a > 0, there is 
a new critical value of F, Fc^ , such that the initial 
datum (C/i(F/yl), 0), which is close to the spiral point 
{zi{F/A), 0), may evolve to the other stable spiral point, 
(z3(F/A),0), for F e {Fc,_,Fc,); see Fig. [H The tra- 
jectory leaving a neighborhood of {zi{F/A),0) and en- 
tering {z^{F/A),0) defines the wave front profile. For 
F < Fc2, the initial datum {Ui(F/A),0) evolves toward 
{zi{F/A),0) and no wave front is generated. Thus Fc2 
yields an approximation to the dynamical Peierls stress, 
Fed- Fig. 1191 compares Fc2 to Fed, which has been cal- 
culated numerically by solving the complete system j^J. 
Notice that our approximation worsens as a decreases to- 
wards zero indicating break down of the one active point 
approximation. We shall explain below why this is so. 

The dynamic critical Peierls stress can be intuitively 
explained as follows. The potential energy associated to 
the nonlinearity h{u) — 2uq + Ag{uo) — F ~ Ui{F/A) — 
UsiF/A) on the right side of Eq. (jTTI) is depicted in Fig. 
I20f a). An initial condition {Ui{F/A),0) is close to the 
left minimum of W{uq). If the energy corresponding to 
the initial condition is slightly higher than that of the 
left minimum, the solution of 117|) evolves toward it be- 
cause of friction. On the other hand, larger energies can- 
not be dissipated by the friction term and the trajectory 
surpasses the maximum of W{uq) with nonzero veloc- 
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ity. Then the trajectory falls in the basin of the right 
minimum of the potential, performing a damped oscil- 
lation about it before reaching Z3{F/A),0). Once uo{t) 
has reached a neighborhood of the second spiral point, 
M_i(t) takes its place and performes a similar motion. 
The resulting trajectory is depicted in Fig. I20f b). The 
solution of the complete system (jJI is shown in Fig. 1^ 
The wave front velocity is approximately given by the 
reciprocal of the time uo(i) takes to jump from zi(F/A) 
to Z3{F/A). Notice that the oscillations about Z3{F/A) 
persist a longer time as a decreases and may have fi- 
nite amplitude. Then we cannot approximate sufficiently 
well ui{t) by the constant value U3{F/A) for small val- 
ues of the friction and the one active point approximation 
breaks down. This explains the discrepancies in Fig. 1191 
for small a. Since the difference [UsiF/A) - Ui{F/A)] 
is larger for the FK than for the cubic nonlinearity, the 
approximation ui{t) « {/^{F/A) is better for the FK non- 
linearity as shown in Fig. 1191 

IV. DISCUSSION 

We have studied wave front solutions in chains of non- 
linear oscillators with inertia and damping. Two ana- 
lytical methods have been used to construct the wave 
fronts and their velocities. For piecewise linear mod- 
els, exact formulas can be found for the wave fronts and 
the relation between their velocity and the applied stress 
F as Atkinson and Cabrera did already in 1965 for the 
conservative case Different from these authors, we 
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FIG. 20: (a) Potential energy W{u) corresponding to the 
cubic nonlinearity h{u) = 2uo + Ag{uo) — F — Ui{F/A) — 
UaiF/A), with y4 = 10 and F = 2.45. (b) Trajectories of Eq. 
llTTt with initial condition {Ui(F/A),0), joining {zi{F/A),0) 
and {z3{F/A), 0) for A = 10, a = 0.5 and F = 2.45. 
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FIG. 21: Numerical solutions Un{t) of Eq. Q with the cubic 
nonlinearity for A = 10, a = 0.5 and F — 5.9. 

have also studied the cases of finite and infinite damping. 
We have depicted the resulting wave front profiles for all 
damping values and found that they may have oscillatory 
tails. For zero damping, these tails oscillate with non- 



decaying amplitude as n — > oo, which means that a new 
definition of wave front is needed 2^. We have shown 
numerically that non-monotone wave fronts with oscilla- 
tory tails (wavy wave fronts) may be stable for certain 
intervals of applied stresses. We have found stable mov- 
ing wave fronts for |F| > Fed > (the dynamic Peierls 
stress), and these fronts coexist with stable static wave 
fronts for an interval Fed < l-^l < Fes {Fes is the static 
Peierls stress; static wave fronts exist for < |i^| < Fes). 

We have also conjectured that the global bifurcation 
diagram for wave front depinning in the presence of iner- 
tia and damping is generically as in Fig. I12r b'). Then 
there are infinitely many saddle-node bifurcations be- 
tween wave front branches in the interval Fed < l-^l < 
-Fes, accumulating at c = and F — Fes, with Fed > 0. 
The basis for this conjecture are our results for the piece- 
wise linear model, that we think are generic also for mod- 
els with smooth nonlinearities given our numerical results 
for them. The function F(c) in the piecewise linear model 
has infinitely many vertical asymptotes (resonances) for 
zero damping that accumulate as c — > and (we sup- 
pose) F — !■ Fes- These resonances become local maxima 
of F{c) as a small damping is added to the model. Max- 
ima and minima of F(c) are limit points (saddle-node 
bifurcations) between branches of wave front solutions in 
the diagram of Fig. ICT b') . 



n 



In the literature |l9l |21 1 , wave fronts with monotone 
profiles have been constructed for undamped models. For 

n 

example, Schmidt [19j showed that w{x) = tanha;, x = 



n + i/2, is a wave front solution for Eq. |5J with a = 
F = 0, A = 1 and a potential V{u) = -Su'^/A - 1% - 
(sinh 1)^2 1^(1^ _ y2 with g{u) = V'{u). For this 

model, i^cs > I . Solving numerically this model for 
F 7^ 0, wave fronts with one oscillatory tail on a; > 
are found for F > 0, similar to those in Fig. ^b). For 
F < 0, wave fronts with one oscillatory tail on a; < are 
found instead; see Fig. ^c). Both types of fronts have 
|c| > 1/2. This shows that wavy wave fronts are generic 
and that Schmidt's monotone wave front is a non-generic 
limiting case separating branches of wavy wave fronts and 
corresponding to Fed = 0. A possible bifurcation diagram 
|c| vs. F for this example could be that in Fig. I12r b) with 
Fes > 0, Fed — and \cm\ — 1/2: numerical solution 
of the model seems to be consistent with a bifurcation 
diagram having a limit point at F = 0, |c| — 1/2, and 
with reflection symmetry with respect to the |c| axis (|c| 
in the diagram is even in F). 
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APPENDIX A: NONLINEAR MODEL HAVING A 
EXPLICIT WAVY WAVE FRONT PROFILE 

It is fairly easy to find smooth g(u) having wave front 
solutions with undamped oscillatory tails provided their 
dynamics is conservative. Let us choose a profile of the 
form: 



-l + fcie*^^ a;<0 
1 + fc2 cos(6x -t- c) a; > 
This profile is continuous if — 1 + fci = 1 + A:2 cos c, 
differentiable if fcio = — /c2^sind and twice differen- 
tiable if k\c? = — A;2&^cosc. We set c = tan^^(6/a), 
fci 2/(1 + a^/62) and fc2 = -A:ia/(&sinc). The only 
restriction on a and 6 is & + tan^^(6/a) < 7r/2, a.h > 
to ensure that cos(6a;-l-c) is monotone in < a; < 1. The 
profile w{x) satisfies; 

a^{w + 1) for w < ki — 1, 

—iP'itv + 1) for w > ki — 1, 
w{x -f 1) — 2w{x) + w{x — 1) = f{w{x)), 

with: 

2(ii; -I- l)(cosha — 1) w<w\ 

(u; -I- l)e"° + y - 2u; (wi,i«2) 

fca cos(VF + 6) + /ciet('^-''-'=) - 2w (w2, t^s) 

2(w — l)(cos5 — 1) w>w-i 

Here, W — k^^ co^^^{w — 1), y = fc2 cosjfoa^-'^ ln[(w -f 
l)//ci] + 5 + c}, w\ = e^°'k\ — 1, W2 = fci — 1, and — 
/c2 cos(5 + c) — 1. The continuous fmiction ](w) is convex 
but, for a certain interval of values of c^, the function 
(7c (w) = c^Uxx — /(w) has three zeroes and is strongly 
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asymmetric. For these values of c, the system (Pun/dt^ = Un{t) = win — ct). 
Un+i — 'iun + Un-i — gc{un) has a traveUng wave solution 
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